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Introduction  and  Review  of  Earlier  Work. 


This  paper  deals  with  estimation  problems  in  which  one  or  more 
observations  are  made  on  a  p-component  vector  X  with  mean  vector 
^  X  =  14  and  covariance  matrix 

(1.1)  '^(X)  =  ^(X-lJ)(X-u)‘ =  I  , 

where  y  and  S  may  have  linear  structure.  The  mean  y  may  be  a 
linear  combination 

r 

(1.2)  u  -  6,  z, 

-  Jit  J  -1 

of  known  p-component  vectors,  z. ,  ...  ,  z  .which  are  assumed  (for 

~i  ~r 

convenience)  to  be  linearly  independent.  The  covariance  matrix 
may  be  a  linear  combination 

m 

(1.3)  Z  =  y  a  G 

g=0  ® 

of  known  symmetric  p  x  p  matrices  Gq,  G^^,  ...  »  »  which  are 

assumed  to  be  linearly  independent;  it  is  also  assumed  that  there  is 

at  least  one  set  ...  »  such  that  (1.3)  is  positive 

definite.  The  coefficients  3,,  ...  ,  3  and  a.,  ...  ,  a  are 

1  r  0  1  m 

parameters . 

If  Z  is  known  or  known  to  within  a  constant  of  proportionality 
and  one  observation  x  is  made  on  X  ,  the  model  is  the  familiar  one 
of  regression  analysis.  The  best  linear  unbiased  estimates  or  Markov 
estimates  of  3j^>  ...  ,  3^  are  the  solutions  to  the  normal  equations 


(1.4) 


X 


9 


r  . 


j=i. 


If  X  has  a  normal  distribution,  (1.4)  are  the  likelihood  equations, 
obtained  by  setting  equal  to  0  the  derivatives  of  the  likelihood 
function  with  respect  to  3^,  ...  >3^  >  and  the  solution  constitutes 
the  maximum  likelihood  estimates.  In  ai;iy  case  the  estimates  are 
unbiased,  ^  3^  =  3^  ,  i=l,  ...  ,  r  ,  and  the  covariance  matrix  of 
the  estimates  is 


(1.5) 


e^)]  -  [z[  I  ^  z.]  ^ 


If  there  are  N  observations  on  X  ,  say  x^,  ...  ,  x^^  ,  the  best 
linear  unbiased  estimates  and  maximum  likelihood  estimates  under  nor¬ 
mality  are  the  solution  to  (1.4)  with  x  replaced  by  the  sample  mean 


(1.6) 


and  the  covariance  matrix  of  the  estimates  is  1/N  times  (1.5). 

Estimation  of  a„,  a. ,  ...  ,  a  was  considered  by  T.  W.  Anderson* 

0  1m 

(1969),  (1970)  when  several  observations  were  made  on  X  and  y  was 
completely  unspecified.  In  this  present  review  we  shall  assume  initially  that 
y  is  known  and  suitably  modify  the  statements  of  the  earlier  papers . 

We  assume  that  X  has  a  normal  distribution  and  that  there  are  N 

observations  x^ . x^^  on  X  .  (N  is  not  necessarily  as  large 

as  p  ;  in  fact,  N  may  be  1  in  some  cases.)  Let 


*  The  1970  paper  was  written  first,  but  there  was  a  delay  of  four  years 
between  its  receipt  by  the  editors  and  its  publication. 


2 


N 


(1.7) 


S  -  if  ' 

a=l 


Then  maximum  likelihood  estimates  of  0.,  a, ,  ...  ,0  are  a  solntion 

0  1  m 

of  the  likelihood  equations 


m 


(1.8) 


^G  =tr( 

m 

iz 

d.  G. 

(2= 

g  , 

~g  \ 

li=0 

h  ~ 

\  h=0 

h 

~hj  ~ 

g=o. 

1. 

•  •  •  ^  HI  y 

these  equations  result  from  setting  equal  to  0  the  derivatives  of  the 

likelihood  function  with  respect  to  0rt,  0, .  0  .  There  is  at 

0  1  m 

least  one  solution  ^1*  ***  *  ^m  (1.8)  such  that 


(1.9) 


m 


S  =  > 


g=0 


d  G 
g 


is  positive  definite.  (The  argviment  given  by  T.  W.  Anderson  in  (1970) 
was  stated  for  C  positive  definite,  but  that  assumption  is  not  needed 
in  general .  If  there  is  more  than  one  solution  to  the 
likelihood  equations,  the  absolute  maximum  to  the  likelihood  function 

A 

is  attained  by  the  solution  minimizing  |2|  .  The  estimates 

Oqs  ^^5  •••  »  consistent  and  asymptotically  efficient  as  ; 

(0^-0^),  (0,-0^),  . . .  ,  (0  -0  )  have  a  limiting  normal  distri 

u  u  1  i  mm 

bution  with  means  0  and  covariance  matrix 


(1.10)  ih  tr  G,,  E  ^  G  . 

~  ~h  ~  ~g 

These  results  follow  from  the  usual  asymptotic  theory  of  maximvim  like¬ 
lihood  estimates. 
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Now  let  us  consider  the  estimation  problem  when  both  y  given 

by  (1.2)  and  Z  given  by  (1.3)  are  unknown.  We  are  to  estimate 

6,,  ....  3  and  a_,  a,,  ...»  a  .  When  X.  is  normally  distri- 

1’  r  0  1  m 

buted  and  ...  »  are  observed,  the  likelihood  equations  are 


(1.11)  511 


i=l 


-1  -1  — 
z.  3.  =  z’.  Z  X 


j=l,  ...  ,  r  , 


(1.12)  tr 


g— 0 ,  1,  ...  ,m, 


where 

(1.13) 


G  = 


N 


a=l 


(1.14)  fi  '  f:  6  ; 

1=1 


and  Z  is  given  by  (1.9).  Then  the  estimates  are  consistent  and 

/V  /s 

asymptotically  efficient  as  .  Moreover,  /n  (3^-3^) »  ...  ,  ^^r”^r^ 

and  Vn  (Oq-Oq)  ,  ’  * "  »  ^^m~^m^  have  a  limiting 

normal  distribution  in  which  the  two  sets  are  independent  and  each  set 

has  the  covariance  matrix  given  previously,  (1.5)  and  (1.10).  CWe 

^  r — ■  N 

note  in  passing  that  asymptotically  as  N-x»  ,  C,  C,  and  (l/N)  ^ 

•V  ^  Oir“  JL 

(x  -x)Cx  -x) '  are  equivalent;  replacement  of  C  in  (1.12)  by  the 
last  matrix  above  represents  a  simplification  of  the  equations.) 


The  relation  between  the  estimation  of  3^^ . 3^  and 

Oq,  a^,  ...  ,  separately  was  indicated  by  T.  W.  Anderson  (1969). 


4 


Consider  (1.4)  for  S  known  and  consider  (1.8)  for  U  known.  In  the 
case  of  normality  the  covariance  between  the  i,j-th  and  k,Jl— th  elements 
of  C  is 


(1.15) 


Gov(c^^,  Cj^g^)  (a^j^  a^j^) 


Let 


(1.16)  c  = 


'11 


'22 


PP 


'12 


,  O  = 


11 


22 


PP 


12 


’  ih  “ 


,(h) 

’ll 

,(h) 

’22 


g, 


(h) 

PP 

(h) 

12 


'P-1> 


P-1.  P 


,(k) 

’p-1,  p/ 


where  G,  = 


h 


^Ij 


Then  c  =  cf  and  (1.3)  can  be  written  as 


m 


(1.17) 


a  = 


h=0 


\  §h  ’ 


which  is  of  the  form  (1.2)  with  U  replaced  by  O  ,  . 3^ 

replaced  by  0^,  0^,  ...  .  Iq’  §1 . §n  ’ 

and  =  y  by  ^c  =  0  .  Then  (1.15)  can  be  written  as 

(1.18)  t  > 

where  $  =  (^^.  j^^)  for  1^3  .  k  _<  £.  and 


5 


(1.19) 


9 


ij,U 


1^ 

N 


i  <  j ,  k  <  Jl  ; 


(1.18)  is  of  the  form  (1.1)  with  X  replaced  by  c  ,  Ji  replaced  by 
a  ,  and  E  replaced  by  #  . 

It  was  shown  by  T.  W.  Anderson  (1969)  that' 


(1.20) 


Ih 


.-1 


c  =  %  tr  Z 


-1 


eh? 


-1 


that  is,  the  bilinear  form  on  the  left-hand  side  of  (1.20)  is  algebrai¬ 
cally  identical  to  the  right-hand  side,  which  is  a  form  appearing  in 
the  likelihood  equations  for  ...  »  *  Substitution  of  g^ 

for  .  c  in  (1.20)  yields 

(1.21)  g;  g„  »  %  tr  Z"^  G,  Z"^  . 

~h  ~  ~g  ~  ~h  -v  ~g 


Thus  the  "normal  equations" 


m 


(1.22)  YZ.  <  8f  =  si  c  , 


=  o, «  ft-l 


h=0,  1,  ...  ,  m  , 


are  identical  to  the  equations 


m 


(1.23)  e.  tr  z“^  G  z”^  G 

i  -V  -vg  -vX 


f«0 


=  tr  Z“^  G  Z“^  C  , 
~  ~g  -  - 


g— 0,  1,  ...  ,m. 


If  the  equations  (1.23)  were  available,  they  would  give  "estimates" 

which  were  linear  in  C  and  unbiased,  and  among  such  "estimates"  they 

would  have  minimum  variance.  (Since  C  is  sufficient  for  Z  ,  these 

would  be  minimum  variance  unbiased  "estimates"  of  a.,  a,,  ...»  a  .) 

0  1  m 
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The  covariance  matrix  of  the  "estimates"  would  be  1/N  times  (1.10). 
However,  these  "estimates"  are  unavailable  since  E  is  unknown;  in 
fact,  the  problem  is  to  estimate  S  . 

The  likelihood  equations  (1.8)  for  . 0^  can  be 

written 

(1.24)  YU  tr  E"-'  G  E"-"  G  8  =  tr  E"-"  G  E"-"  C  ,  g=0,  1,  ...  ,  m  , 

-  ~g  ~  I  ~  ~g  ^  ~ 

^  1 

by  multiplying  the  left-hand  side  of  (1.8)  by  I  =  V  8^  G^  E 

These  .equations  are  similar  to  (1.23),  but  E  in  (1.23)  has  been 
replaced  by  E  .  As  will  be  shown  later,  the  form  (1.23)  suggests 
computational  procedures  and  asymptotic  properties.  ■ 

One  of  the  probability  models  in  which  the  covariance  matrix  has 
the  form  (1.3)  is  a  moving  average  stationary  stochastic  process  of 
finite  order.  Let 


m 


(1.25) 

X  =  a 

,  ‘  ® 

Vg  » 

,  1,  0,  1,  . . ,  , 

where 

=  0  ,  ^v^  =  0^  , 

and  $v^  V  =  0  ,  t^s 
^  t  s  ’ 

.  Then  x^  =  0 

and 

(1.26) 

^h 

-h  -  ^+h 

2 

“  ^  V  ot  a  , ,  , 

^  s  g-H.  ’ 

h— 0,  1,  ...  ,  m  , 

II 

o 

h  =  m+1 . 

The  vector  x  =  (x  ,  ...  ,  x  )’  has  the  covariance  matrix  (1.3)  with 
"  .  1  P 


7 


•1 


When  cTq,  a^,  ...  ,  are  defined  by  (1.26)  real  aQ,  ...  .  > 

2 

^  0  >  2  given  by  (1.19)  is  positive  definite  (of  any  order). 


8 


9  ®  •  •  9 


If  X  =  (x^ 

of  ^1*  • • •  > 

(for  y  =  0  and 

is 


x  ) ’  is  observed,  the  maximum  likelihood  estimates 
P 

0^  are  defined  by  (1.8)  or  (1.24)  where  C  =  xx' 

N  =  1).  Then  the  right-hand  side  of  (1.8)  and  (1.24) 


(1.30) 


-1 

tr  E  G 
~g 


x’  E  G  E  X  . 

-  "  -g  ~  - 


In  the  case  of  the  moving  average  model  one  may  be  interested  in 

2 

the  parameters  a^,  ...  ,  ,  0  ,  with  =  1  .  If  the  estimates 

»  •••  >  S  3.ve  such  that  the  estimated  covariance  matrix  obtained 
0  1  m 

by  replacing  0q,  0^^^,  ...  ,  0^  in  (1.29)  by  . . .  ,  is 

positive  definite  for  (1.29)  of  every  order,  then  the  estimate  of  the 
spectral  density 

1 

(1.31)  f(X)  -  y  ^  0  cos  Xg 

is  positive  and  the  equations  (1.26)  with  0q,  0^,  ...  »  replaced 

2 

by  8_,  8,,  ...  ,8  and  a,,  ...  ,  a  and  0  replaced  by  8, ,  . . .  ,  8  , 

■^0  1  m  1  m  ^  ■'I  m 

2  2 
8  can  be  solved  for  8^,  ...  ,  8^  ,  8  >  0  .  The  polynomial  equation 

associated  with  the  moving  average  (1.25)  is 

m 

(1.32)  ni  a  =  0  . 

g=0  s 

If  the  roots  are  required  to  be  not  greater  than  one  in  absolute  iralue, 

2 

a. . a  and  0  are  uniquely  determined  by  0^,  0, ,  ...  ,0 

1  m  0  1m 

[See  Section  5.7  of  T.  W.  Anderson  (1971). J 
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In  this  case  of  the  time  series  problem  where  N  =  1  ,  we  are 
interested  in  the  as37mptotic  theory  when  p-^  .  These  properties  will 
be  studied  later. 

It  was  pointed  out  by  T.  W.  Anderson  (1969),  (1970)  that  the  model 
described  above  is  appropriate  for  many  problems  of  the  analysis  of 
variance.  For  example,  let 

(1.33)  x^^  =  y  +  (yj.-y)  +  u^  +  •••.?» 

0t“l  f  y  ^  y 

where  y  =  }  y  /p  ,  ^u  =  0  ,  '^v,  =  0  ,  %,u  =  a  ,  :g,v  =  a  » 

and  all  “qj's  and  independent.  Then,  ~  y^  the 

covariance  matrix  of  x,,  ,  . . .  ,  x  is  given  by  (1.3)  with  =  I 

xOt  pCi  -- 

and  G-  =  ee'  ,  where  e  =  (1,  ...  ,1)*  »  This  is  a  mixed  model  in 

the  analysis  of  variance.  The  overall  mean  is  y  ,  the  (y^-y)'s  are 

the  fixed  factor  effects  and  the  u  's  are  the  random  factor  effects. 

a 

The  factor  analysis  model  when  the  factor  loadings  are  knowm  was 
studied  by  T.  W.  Anderson  (1970).  This  model  is  particularly  appro¬ 
priate  for  one  form  of  Guttman's  quasi-simplex. 

Hartley  and  Rao  (1967)  have  given  the  derivative  equations  (1.4) 
and  (1.8)  for  N  =  1  when  G«,  G, ,  ...  ,  G  are  generated  by  models 
of  the  analysis  of  variance.  Their  proposals  for  solution  are  different 
from  the  one  presented  by  T.  W.  Anderson  (1970)  and  the  one  presented  in 
Section  2  of  this. paper.  Other  references  were  given  by  T.  W.  Anderson 
(1969),  (1970). 
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2. 


d 


_ Computation  of  Estimates  for  Covariance  Matrices 

The  set  of  equations  (1.24)  suggests  an  iterative  method  of  solving 
the  likelihood  equations  for  0^,  0^,  ...  ,  .  Let  •••  » 

be  an  initial,  set  of  values;  these  may  be  values  given  a  priori  or  they 
may  be  estimates  obtained  in  another  way.  Let  •••  > 

the  solutions  to 


(0) 


m 


(2.1) 

where 

(2.2) 


f=0 


tr  sT\  G  G  d  =  tr  2“^  G  E  C  ,  g=0,  1, 

~i-l  -g  ~i-l  ~f  f  -i~l  -g  ~i-l  ~  ^  » 


•  •  •  «  IQ 


m 


1  -  - 


i=l.  2. 


i=l,  2. 


•  •  •  t 


The  equations  (2.1)  can  also  be  written 
(2.3) 


m 


i;  il-i  If  »£  - Ji-i  s  ■ 


g=0>  1 . . 


where  ,  is  formed  from  E.  ,  as  <E>  is  formed  from  E  . 
~i-^l  ~i-l 


m 


Lemma  2.1.  If  2 _ *^h  ^h  is  nonsingular 


(2.4) 


tr 


m 


h=0 


“h  ?h)'^  2.  (£  2.)-^  2f] 


is  positive  definite. 

Proof .  For  (y^,  y^,  ...  ,  y^)  (0,  0,  ...  ,  0) 


m 


(2.5)  tr  _  a,  G  f’5  a,  G.V^  G 

gTf^o  vfco  ^  ^  ~ 


'=£  ='g  ^£ 


tr 


Yf  9 


>  0 
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9 


because  Gq,  G^,  ... 


G  are  linearly  independent.  Q.E.D. 


-m 


If  ^  is  positive  definite,  the  matrix  of  the  coefficients  of 
Qqj  ^2.’  * '  ”  ’  ^m  the  left-hand  side  of  (2.1)  is  positive  definite, 
and  hence  there  is  a  unique  solution.  (The  iterative  procedure  suggested 
here  can  be  considered  as  an  approximation  to  the  procedure  proposed  by 
T.  W.  Anderson  (1970)  on  page  6;  the  present  proposal  is  computationally 
simpler  and  its  properties  can  be  studied  more  easily.)  The  iteration 
may  be  stopped  at  the  i-th  stage  if  ^ .  m  ^  does  not 


differ  by  much  from  : 


(i-l) 


m 


Since  '^C  =  2  ,  given  by  (1.3),  unbiased  estimates  of  Oq,  O^, 


O  can  be  obtained  as  the  solutions  to 
m 


m 


(2.6) 


f=0 


S.  tr  0  G  0  G-  =  tr  0  G  0  C  , 

f  - g  ~  ~f  ~  ~g - 


g— 0 ,1,  ...  ,m. 


for  an  arbitrary  positive  definite  matrix  0  .  These  estimates  (under 
normality)  have  covariances 


(2.7)  t.0,.  i 


m 


N 


h,j ,k,£=0 


m^^  m®^  a,  a.  tr  0  G,  0  G.  0  G.  0  G,  , 
k  a  ~  ~h  ~  ~J6  ~  ~3  ~  ~k 


f ,  g— 0 ,1,  ...  ,m. 


where  (m^^)  =  (m„  )  ^  and 
rh 


(2.8) 


®fh  =  ®  ®  5h  ’ 


f ,  h=0 ,  1 ,  ...  ,  m  . 


As  these  estimates  are  consistent  and  ^  ’ 

(0  -0  )  have  a  limiting  normal  distribution.  [If  Ji  is  unknown,  C 
mm  s.  ’  ~ 

in  (2.6)  could  be  replaced  by  [N/(N-1)]C.] 
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The  equations  (2.6)  can  be  obtained  from  (1.23)  by  replacing  S 
by  0  .  A  particular  choice  of  0  is  I  ;  this  substitution  corres¬ 
ponds  to  the  Markov  estimates  when  the  variances  of  /N  c^^  are 
proportional  to  2,,  the  variances  of  ^  ^  J  proportional 

to  1  and  every  covariance  is  0  . 

To  obtain  asymptotically  efficient  estimates  of  *^1’  *  *  *  ’  '^m 

only  one  step  in  the  iteration  is  needed  if  the  initial  estimates  are 

consistent..  See  Section  4.  If  N  >  1  and  y  is  unknown,  C  may  be 

- - -  _  _ 

replaced  by  (1/N)  } _ m-i  the  computation;  one  may 

wish  to  multiply  the  solution  by  the  factor  N/ (N-1)  . 

The  solution  of  (2.1)  requires  evaluation  of  quantities  such  as 

(2.9)  tr  A“^  B  A“^  L  , 

where  A,  B,  and  L  are  symmetric  and  A  is  positive  definite,  finding 
A  ^  corresponds  to  solving 

(2.10)  AX  =  I  . 

The  "forward  solution"  of  a  method  of  pivotal  condensation  or  successive 
elimination  corresponds  to  multiplying  (2,10)  on  the  left  by  a  triangular 
matrix  F  to  obtain 

(2.11)  TX  =  F  , 
where 
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(2.12) 


(2.13) 

Then 

(2.14) 

that  TF ' 
T  and  F 
above  the 
symmetric 
Let 


0 

0 


“  \ 


0 


has  O’s  below  the  main  diagonal  follows  from  the  facts  that 
have  O’s  below  the  main  diagonals,  and  that  TF’  has  O’s 
main  diagonal  follow  from  the  fact  that  TF’  =  FAF’  is 

Since  FAF’  is  positive  definite,  t..  >  0,  i=l,  ...  ,  p  . 


14 


/ 


(2.15)  D  = 


/t 

0 

0 


11 


0 


22 


0 

0 


33 


and  let  H  =  D 

(2.16) 

(2,17) 


Then 


0 

0 

0 


/r~  / 

pp  / 


HM’  =  I  , 
a“^  =  H’H  . 


(This  development  is  given  in  more  detail  by  T*  W.  Anderson  C1971)  in 
Section  2.3.)  Note  that  only  the  forward  solution  is  needed  to  obtain 

H  . 

Then 


(2.18)  tr  A“^  B  a“'^  L  =  tr  H'H  B  H'H  L 

=  tr  H  B  H'H  L  H'  . 

Thus  the  symmetric  matrices  HBH'  and  HLH'  are  computed.  When  N  =  1 
and  y  =  0  ,  C  =  xx'  and  (2.18)  with  L  replaced  by  C  =  xx'  becomes 

(2.19)  tr  H'H  B  H'H  xx'  =  x'  H'H  B  H'H  x 

=  (Hx)  'HBH'  (Hx) 

=  (Fx)'  F  B  F'  CEx) 

=  (F'  Fx)  B  (F'  d“^  Fx)  . 

/S 

Given  A  =  ,  the  data  enter  the  equations  (2.2)  through  HGH*  . 

In  the  case  of  C  =  xx'  ,  this  involves  only 
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(2,20) 


9 


Hx  =  d"^  F  X 


which  is  the  forward  solution  applied  to  x,  followed  by  the  division 
of  each  element  of  the  resulting  vector  by  the  square  root  of  the 
diagonal  element  of  FA  =  T  . 

If  (1.2)  holds  and  •••  ,  3^  are  unknown,  initial  unbiased 

estimates  of  3- ,  ...  ,  3  can  be  obtained  from  the  equations 
X  IT 

r  /V  _ 

(2.21)  YU  2!  0  Z  3.  =  z!  0  X  ,  j=l . r  , 

i=l  ~  i  ~  ~ 


where  0  is  any  positive  definite  matrix.  The  solution  has  a  multi- 
variate  normal  distribution  with  covariance  matrix 


(2.22) 


i  (2*  0  z.)"^(z’  0  S  0  2.) (2’  0  2,)"^ 


These  estimates  are  consistent  as  .  If  0  =  I  ,  these  estimates 

are  least  squares. 

When  32^»  •••  >  3^.  and  Oq,  ...  »  are  unknown,  the  esti¬ 

mates  obtained  from  (2.21)  can  be  inserted  into  (1.5)  to  obtain  an 
estimate  of  y  and  this  estimate  in  turn  can  be  used  to  define  G  . 

Then  O. ,  ...  ,  O  can  be  estimated  and  E  :  this  estimate  of 

E  can  replace  S  in  (1.4)  to  obtain  improved  estimates  of  3^ .  3^ 

This  procedure  yields  consistent  and  as3miptotically  efficient  estimates 

of  3 . ^  j  •••  >  ^  ns  .  We  shall  study  later  their 

1  r  O  l  m  ^ 

asymptotic  properties  as  p-^  . 
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3»  Computation  for  the  Moving  Average  Process 

3.1.  Case  of  First-Order  Moving  Average  Process.  The  special 
form  of  E  as  given  in  (1.26)  makes  the  computation  easier  in  the 
case  of  estimating  the  nonzero  covariances  of  a  finite  moving  average 
process.  This  is  of  particular  interest  in  the  case  of  one  observa¬ 
tion  X  when  y  is  assumed  0  .  To  illustrate  the  procedure  we 
consider  the  case  of  order  m  =  1  .  Let  p  =  »  and  let  the  p  x  p 

matrix  (1/CT^)E  be 


(3.1) 


All  of  the  elements  of  this  matrix  are  0*s  except  on  the  main  diagonal 
and  one  above  and  below  the  main  diagonal.  The  method  of  pivotal  con¬ 
densation  or  successive  elimination  starts  with  leaving  the  first  row 
unchanged  and  subtracting  p  times  the  first  row  from  the  second  row; 
this  operation  changes  the  second  row  to  a  row  having  nonzero  elements 
only  on  the  diagonal  and  one  entry  to  the  right  of  the  diagonal.  Each 
successive  step  consists  of  subtracting  a  suitable  multiple  of  one  row 
from  the  next.  This  "forward  solution"  can  be  represented  as  FA, 
where 


(3.2) 


F  =  F  F  ,  .. 
~  -p  ~p-l 


F  F 
-3  ~2 


(3.3) 
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The  matrix  is  the  identity  except  that  the  element  in  the  j-th 

row  and  (j-l)-st  column  is  -p  times  the  reciprocal  of  the  element 

in  the  (j-l)-st  row  and  (j-l)-st  column  of  F .  ^  F .  »  . . .  F,  A  . 
Thus  F  has  the  form 


where  the  orders  of  the  I*s  are  j-2  and  p-j  ,  respectively.  (In 
the  product  F  =  F^  ^p-1  * ’ '  ^2  ^1  element  j  i  appears  in 

the  j-th  row  and  (j-l)“St  column.)  The  calculation  of  Hx  involves 
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(3.7) 


9 


,  Fx  =  F  F  T  . . .  F_  F-  F.  X  =  w 
~,p  ~p-l  ~3  ~2  ~1  ~ 


say.  The  computation  of  w  proceeds  as  follows:  » 

(3.8)  “J  -=‘j  •  J-2. 


Thus  the  elements  of  w  can  be  calculated  in  sequence . 

33 


We  can  write  an  equation  for  f,  .  .  Let  be  the  j,j-th 

J » J 


term  of  'F,  F,  ,  .  . .  F.,  A  .  The  method  of  successive  elmination 
~3  ~1  ~ 

shows  that 


(3.9) 


(3.10) 


j+l,j  ’ 

jj 

a(j+l)  =  1  _  _pL  _±i _ 1. 

^j+l,j+l  ^  (j)  (j)  • 

jj  jj 


The  elements 

in  sequence. 

Next  let 

(3.11) 

-2 

D  w  =  u  , 

where  D  is 

defined  in  (2.15);  in  components  this  is 

(3.12) 

u  =  w,/d^,  =  w  /t..  =  w,  /af^^  , 
j  J  JJ  j  JJ  J  JJ 

Finally,  let 

(3.13) 

V  =  F'  u  =  F'  F'  ..  F'  F'  w  ; 

~2  ~p— 1  ~p  ~  ’ 

thus  v'  =  u 
•  P  P 

1 

(3.14) 

V,  =  u,  +  f.  ...  .  V.  ,  j  =  p-1, 

J  J  J+1>J  J+1  ’  ’ 

j— 1,  ...  ,  p 
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Then 


(3.15) 


i'  ?h  i  ^  5h  ? 

“o 


■  3  x'  Sh  X 


h=0,  1  . 


For  =  I  ,  V  G»  V  =  v*v  =  )  v  ,  and  for  G  given  by  (1.27) 

-V  -vU  ^  -w  -w  .  3“-^  J  ' 

(3.16) 


y’  2i  y "  ^  ^j+i  • 


The  number  of  arithmetic  operations  in  the  calculation  of  these  quad¬ 
ratic  forms  is  approximately  proportional  to  p  , 

To  calculate  the  coefficients  of  the  unknowns  in  the  iterative 
procedure  we  need 


(3.17) 


tr  A“^  G  A“^  G,  =  tr  H’H  G  H’H  G, 
~g  ~  ~h  - g - h 


=  tr  F’  D~^  F  G  F’  D~^  F  G, 

-  -  -  ~g  -  ~ 

=  tr(D'“^  F  G  F’  D"^)(d"^  F  G,  F’  d“^) 


The  forward  procedure  F  can  be  applied  to  each  G^  on  the  left  and 

F'  on  the  right  followed  by  D  ^  on  the  right  and  left  to  obtain  the 

—1  —1 

S)raimetric  matrix  D  F  G  F*  D  .  However,  in  this  case  of  m  =  1  , 
-1  -1 

the  form  of  E  =  (l/aQ)A  is  known  and 


(3.18) 


tr  G  Z"^  G,  =  -4-  tr  A“^  G  A“^  G, 
~  ~g  -v  ~h  2  ~  ..h 

°0 


can  be  computed  directly.  See  Shaman  (1969). 


The  matrix  G^  can  be  written 


(3.19) 


Si  -  r  - 
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where  the  orthogonal  matrix  P  is 
(3.20) 


=  / 


p+1 


)  ■ 


and  the  diagonal  matrix  A  has  2  cos  TTj/(p+l)  as  its  j-th  diagonal 
element.  [See  Section  6.5.4  of  T.  W.  Anderson  (1971).]  Then  the  quan¬ 
tities  (3.18)  can  also  he  written 

(3.21)  tr  f  ^  i  tril  +  p  G^]'^  Sg  Ij  +  P  5h 


a. 


^  tr[PP'  +  P  P  A.  P*]"^  PAP'  [PP'  +  p  P  At  P']“^  P 


^0 
1 


=  tr[P(I  +  p  A^)P’]"^  ?  ^  +  P  "  I  hil' 


-1 


a 


■4  tr  P(I  +  p  A, )"^  P’P  A  P’P(I  +  p  At)“^  P’P  a.  P 
2  ~  ~  ~g  ~  ~  ~  ~  ~  ~n  ~ 


I  tr(I  +  p  Aj^)“^  Ag  (I  +  p  A^)"^  A^  , 


g,  h=0,  1, 


where  Aq  =  I  .  Then 


(3.22) 


(3.23) 


(3.24) 


“  6"^  ?o  ^  So  ■ 


tr  a“^  G-  A“^  Gt  = 


tr  a“^  Gt  a"^  Gt  = 

'  ^  '^»X  '>#X 


j=l  (l+2p  cos  ^  )^ 


i!: 


2  cos 


p+1 


j=l  (l+2p  cos  )^ 


4  cos 


2  T[j_ 
p+1 


j=l  (l+2p  cos  ^  )^ 


These  sums  can  be  approximated  by  integrals.  The  sum  (3.22)  is  approxi¬ 
mated  by 
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i  > 


(3.25) 


£±1 

TT 


IT 

2(p+l) 


ir 

2(p+l) 


dX 

(l+2p  cos  X)^ 


which  in  turn  is  approximated  by 


(3.26) 


2±1 


TT 


/: 


dX 


0  (l+2p  cos  X)^  2(p+l)  L(n.2p)2  (l-2p) 


=  £±1  f 

Jr 


TT 


dX 


1+Ap 


0  (l+2p  cos  X)^  (l-4p^)^ 


2  3/2 

The  first  term  on  the  right-hand  side  of  (3.26)  is  (p+l)/(l-4p  ) 

_2 

[See  Pierce  (1929),  Formulas  300  and  305,  for  example.]  Then  tr  E 
is  approximated  by 


(3.27) 


A 

4 


(p+1)  (I4a^)  _  n-6a^  + 

(1-a^)^  (1-a^)^ 


1  p  -  6a^  -  (p+2)a^ 
(l-a^)^ 


In  a  similar  way  (3.23)  and  (3.24)  are  approximated  by 


(3.28) 


4(p+l)p 

(l-4p^)^/^ 


8P 


(l-4p2)^ 


(3.29) 


4  (p+1) 

2 

P 


1  - 


l-8p 


-  4 


l+4p 


2  2 
(l-4p^)^ 


4p  , 


P  0  , 

p  =  0  . 
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The  2x2  matrix  of  coefficients  tr  S  ^  S  ^  (g,  h=0,  1)  is 


approximated  by 


(3.30) 


p+1 


o'*  (l-a2,3 


/  I-KX^ 


-4a 


-4a 

2  4 

2+8a  -  2a 


l+6a^  +  a^  -8a(l4a^) 


(1-a^)^  I -8a(l+a^)  4(1-Ki^) 


An  alternative  approach,  which  may  be  generalized  to  cases  of  m  >  1  , 
is  to  approximate  the  moving  average  covariance  matrix  by  the  inverse  of 
an  autoregressive  covariance  matrix.  The  covariance  matrix  of  the  moving 
average  process  of  order  1  may  be  written 


/ 


l+a 

a 

0 


0 

0 


14a' 


a 


0 

0 


0 

a 

I4a^ 


0 

0 


0 

0 

0 


14a 

a 


0 

0 

0 


a 


1+a^  j 


The  matrix  in  the  exponent  of  the  first-order  autoregressive  Gaussian 
process  satisfying  the  stochastic  difference  equation 


(3.32) 


yt  +  ^‘yt-i  =  \  . 


^9  ^9  ^9  •••  y 
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IS 


(3.33) 


'i' 

Iar 


u 


a 

r 

l-Kx'' 

a 


0 

0 


0 

a 


0 

0 


0 

0 

0 


1-HX 

a 


0  \ 
0 
0 


a 

1 


This  differs  from  E  . 

-vMA 


2  2 

by  a  a  E  ,  where 


(3.34) 


E  = 


h 
0 

0 


0 

0 


0 

0 

0 


0 

0 


0 

0 

0 


0 

0 


0 

0 

0 


0 

0 


if  =  1/a^ 
u 


This  fact  suggests  approximating  by 


/ 


5ak  =  !S- 


1-a 


1 

-a 

2 

cx  •  •  • 

-a 

1 

•  •  • 

2 

a 

• 

-a 

• 

1 

• 

• 

• 

• 

• 

• 

• 

• 

• 

(-a)^  ^  ... 

• 

• 

1 

• 

• 

-a 

(-a)^~^ 

(-a)^~^ 

(-a)^  ^  . . . 

-  a 

1 

with  =  1/a^  ;  Then 
u 
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1 


(3-36)  =  ^ 


l4a^-a^-a^P 

-a[2-a^-a^P 

a^[3-2a^-a^P“^J 

-a[2-a^-a^P"^] 

14a^-a^-a^P“^ 

-a[2-a'^-a2p-'^J 

a^[3-2a^-a^P"‘^] 

-a[2-a^-a^^ 

_  2  6  2p-4 

14a  ,-a  -a 

• 

, 

• 

\(-a)P"^[p-(p-l)ci^-a^] 


\ 


2  2  2pl 

14a  -a  -a 


and' 


(3.37)  tr  2^  Gq  2^  Gq  - 


a’  (1-a^)^ 


{p(l-a^)-2a^+2a^P'‘'^}  . 


The  other  quantities  needed  2^^^  G^  2^  ,  tr  2^^  G^^  2^  Gq  ,  and 

tr  2._  G,  2._  G-  can  also  be  computed. 

~M  ~1  ~AR  ~1 

^'-l 

In  the  iterative  procedure  is  replaced  by  2^  where 

a^(l4a^)  is  &Q^  and  a^a  is  .  In  situations  where  p 


is  .  large  and  we  consider  limits  as  p-x»  ,  we  can  use 

1 


(3.38) 


(3.39) 


(3.40) 


lim 

p-^oo 


lim 

p-x» 


lim 

p-w 


tr  2.„  G_  2.„  G_ 
~AR  ~0  ~AR  ~0 


P 

4  2.3 

o  (1-a  ) 

rr  2 , 

G 

2, 

G 

....  2  4 

~AR 

~1 

~AR 

~1 

24-8a  -2a 

P 

(1-a^)^ 

tr  2,„ 

G., 

2. 

G 

~AR 

~0 

~AR 

~1 

-4a 

2.3 

a  (1-a  ) 
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9 


If  the  matrix  A  given  by  (3.1)  is  used  in  place  of  E  =  ,  then  a 

can  be  replaced  by  1  and  a/(l-kx^)  by  p  in  (3.38),  (3.39),  and  (3.40) 
to  give  l/(l-4p^)^^^  and  the  coefficients  of  p+1  in  (3.28)  and  (3.29), 
respectively.  These  values  agree  with  the  limit  of  1/p  times  (3.30). 


(3.41) 


where 


T-h 


(3.42) 


S  -  T-h  ^  "'j+h 

x'  G,  X 
_  ~  ~h  ~ 

2 

tr 

~h 


This  is  a  consistent  estimate  of  p  as  p-^  ,  The  estimates  Cq 
and  are  the  unbiased  estimates  of  Oq  and  obtained  from 

(2.6)  with  0=1. 


3.2.  Case  of  Higher-Order  Moving  Average  Process.  Let  p^^  = 
and  let  the  p  x  p  matrix  (l/aQ)E  be 


(3.43)  A  = 


m 


1 

Pi 


m-1 


m 


m-2 


m-1 


m 


m-1 


m-2 


1 

Pi 


m 


P. 


m-1 


0 

0 

0 
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which  has  0*s  m+l  and  more  positions  above  and  below  the  main  diagonal 
(m  <  p)  .  Consider  again  the  forward  solution  of  a  method  of  pivotal 
condensation  or  successive  elimination.  ¥e  represent  this  as  the  multi¬ 
plication  of  A  on  the  left  by  F  =  .  F2  ,  where  F^^  =  I  . 

In  one  procedure  F2  represents  subtracting  times  the  first  row 

of  F-  A  =  A  from  the  second,  P-  times  the  first  row  from  the  third, 
...  ,  and  p  times  the  first  row  from  the  (m+l)st.  Then  F„  F,  A 

has  all  O’s  in  the  first  column  below  the  first  entry.  Before  the 
j-th  step  F^  ^  ...  F^  A  has  all  O’s  in  the  first  j-1  columns 
below  the  main  diagonal.  F^  represents  subtracting  multiples  of  the 

(j-l)st  row  of  F,  .  . . .  F.  A  from  the  j-th  ,  (j+l)st . (j+m-l)st 

— j—i  ~i  ~ 

rows  in  order  that  the  j-th  column  have  only  O's  below  the  main 
diagonal.  Let  a^^  be  the  i,  k-th  element  of  A^^  ”  ^j  1  *"* 

Then 


(3.44) 


F.  = 
~3 


I'l 

0 

0 

0 


0 

1 

f . 

~3 

0 


0 

0 

I 

0 


0  \ 
0 
0 

h 


where  the  I's  are  of  order  j“2,  m,  and  p+l-j-m  ,  respectively,  and 

/  0.  \ 


(3.45) 


f.  =  - 
~3 


,(j-l) 
j-1, j-1 


m 

(j-1) 

^j-l,j-l 


j=2,  ...  ,  p-m+l  ; 
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for  j  -  p-m+2,  . e .  ,  p  ,  f j  consists  of  the  first  p-j+1  components 
of  (3.45).  Thus  F  consists  of  the  identity  with  (at  most)  m  non- 
zero  elements  below  the  main  diagonal  in  the  (j-l)st  column. 

The  operations  can  be  done  in  another  order.  Then  F2  represents 

subtracting  times  the  first  row  of  A  =  A  from  the  second;  F^ 

represents  subtracting  appropriate  multiples  of  the  first  two  rows  of 

^2  ^1  ^  from  the  third  row  to  obtain  O's  to  the  left  of  the  diagonal 

elements  in  the  third  row  of  F^  F2  F^^  A  .  Before  the  j-th  step 

F,  -  . . .  F,  A  has  O's  in  the  first  i-1  rows  to  the  left  of  the  main 

diagonal  and  the  remaining  p-j+1  rows  of  F^  ^^  .  . .  F^^  A  are  the  same 

as  those  of  A  .  F.  represents  subtracting  appropriate  multiples  of 

the  (j-m)th,  (j-m+l)st,  ...  ,  (j-l)st  rows  of  F^  ^  . . .  F^^  A  in  turn 

to  make  the  (j-m)th,  (j-m+l)st,  ...  ,  (j-l)st  elements  of  the  j-th 

row  of  F ,  . . .  F.  A  zero.  Thus  F,  consists  of  the  identity  with  (at 
-J  ~  ~J 

most)  m  nonzero  elements  to  the  left  of  the  main  diagonal  in  the  j-th 


row. 


The  calculation  of  Hx  is  done  by 


(3.46) 


w  =  Fx  =  F  F  ,  . . .  F„  F,  X  . 

~  ~p  ~p-l  ~2  ~1  ~ 


Let 


(3.47) 


w^^^  =  F, 


Then  w  =  w^^^  .  In  the  first  sequence  of  operations  in  (3.47) 


(i-1) 

w^*^  for  i=l,  ...  ,  j-1  , 
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(3.48) 


^i,j-l  ^j-1  * 


l.“j  y  «  •  a  y  j“Hll““l  y 


w 


=  wp  ,  i  =  j4m,  ...  ,  p  ;  if  j  =  p-m+2,  ...  ,  p  (3.48) 


holds  for  i= j ,  . . .  ,  p  .  Thus  w 
Next  let 


(j) 


•jiaputed  successively. 


(3.49) 


-2 

u  =  D  w 


and 


(3.50) 


V  =  F'  u  =  F'  F'  ...  F'  F*  u 

~  ~2  ~p~l  ~p  ~ 


The  operation  (3.50)  can  be  done  sequentially 


(3.51)  » 

where  =  u  and  =  v  .  Then 


j=p,  ...  »  1  » 


(3.52) 


X*  X 


1  t  n 

-tt  X  A  G,  A  X 

"o 

i  r  Sh  ^  • 


b-*0  J  1)  •••  ^lUa 


?0 

etc. , 

(3.53) 


=  I  ,  v'  G„  V  = 


P  .2 


Gq  V  =  I _ and  for  G^^  given  by  (1.27),  (1.28), 


v’  G,  V  =  2  V.  V, , 

-  ~h  ~  3  j+h 


h"! , '  . . .  ,  n  . 


The  number  of  aritiunetic  operations  is  approximately  proportional  to 
mp  . 
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I 


As  indicated  for  the  case  nii=l  the  coefficients  of  8^,  d,,  ^ 

u  1  m 

ran  be  Calculated  from  the  forward  solution  according  to  (3.17)  and  (3.18). 

The  approximations  to  tr  A  ^  A  ^  presented  for  the  case  m=l 
can  be  extended.  The  matrix  G^  is  approximately  2  A^  treated  in 
Section  6.5.4  of  T.  W.  Anderson  (1971).  Then  G^^  is  approximately 
P  P  ,  where  A^  is  a  diagonal  matrix  with  2  cos  hjTT/(p+i)  as  its 
j-th  diagonal  elements.  (P’  P’  =  P’P  =  I  ,  P’  G,  P*  =  A.  diagonal 
and  P'  G,  P  is  A,  diagonal  plus  a  matrix  relatively  small,  j=2,  ...  ,  m  .) 
Then  we  have  the  approximations 


(3.54)  tr  A"^  G-  A"^  G- 


m. 


(1  +  2  )  p,  cos 
h=l 


h  p+1 


p+1 

TT 


ir  - 


ir 


2  (p+1) 


IT 


dX 


m 


2  (p+1)  (1  +  2  cos  Xh) 


(3.55)  tr  A"^  G«  A"^  G 

~  ~  ~s 


cos 


.  2 


EiEL 

T+1 


(1  +  2  fz  cos  ^  )^ 


~  2 


p+1 


IT 


IT  - 


TT 


2 (p+1) 
TT 

2  (p+1) 


cos  Xg  dX 


m 


(1  +  2  IZ  cos  Xh)‘ 
h=l 


g— 1 ,  ...  ,  m  , 
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cos 


(3.56)  tr  a“^  G.  A"^  G 
~  -f  ~  ~g 


m  cos  ^ 
T+1  T+1 


4  ii  — 

A..I  ■  ■  -  m 

(1  +  2  JZ  cos  ^ 
h=l  P+l 


TT  - 


TT 


~  4 


p+1 


TT 


2(p+l) 


TT 


2(P+1) 


cos  Xf  cos  Xg  dA 
m 


(1  +  2  2 _  cos  Ah) 

h=l 


f )  §~1  j  •  •  •  ,  in  . 

In  the  integrals  .  cos  k  can  be  linritten  as  a  polirnomial  in  cos  A  of 
degree  k  , 

As  in  the  case  of  m=l  ,  the  covariance  matrix  2^^  of  the  moving 

average  process  of  order  m  can  be  approximated  by  ,  the  matrix 

in  the  exponent  of  the  normal  distribution  of  the  Gaussian  autoregressive 

-1  -1 

process  of  order  m  .  Then  2^  is  approximated  by  '1'^  =  2^  ,  whose 

elements  are  the  covariances  of  the  autoregressive  process.  If  the 

roots  of  (1.32)  are  different,  say,  z. ,  ...  ,  z  then  the  i,i-th 

1  m 

element  of  2,-  can  be  written 
~AR 


m 


(3.57) 


\  "h 

h=l 


i-J 


for  suitable  constants  k^,  ...  »  •  [See  T.  W,  Anderson  (1971), 


AR 


j=l 


Then  the  i , 

i-th 

2 

element  of  2 . „ 
-AR 

i-1 

m 

k  k,  z^-j  zj-j 

gag  h 

(i-3)  =  H 
3=1 

IZ 

g,h=l 

m 


+  II  m 


4 
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and 


gVi  i-Vh  - 


<3-59)  tr  4  -  t 
i»j=l 


m 


-  W 


g.li=l 


P(i  +  Vh’  ■  ^  “g^h  m-<^gV‘’J 


1  -  z  z, 
g  h 


Then 


(3.60) 


g,h=l 


r,-4oo 


g 


1  +  Z  Z, 
1 

;  n  1  -  z  z, 
g  ^ 


-1  -1 

This  gives  a  method  for  approximating  tr  E  G  E  G. 

~  ~g  ~ 


i^h  • 


If  the  integrals 


in  (3.54),  (3.55)  and  (3.56)  are  over  the  interval  -n,  rr  ,  they  are 

identical  to  the  evaulation  of  tr  E._  G  E.„  G.  for  the  integrand  in 

~AR  ~g  ~AR  ~h 

(3.54)  is  proportional  to  the  square  of  the  spectral  density  of  the 
autoregressive  process.  (See  Theorem  8.3.3  of  T.  W.  Anderson  (1971).] 
We  primarily  want  to  use  this  approximation  to  obtain  the  coeffi¬ 
cients  of  •••  $  (2.1)  on  the  basis  of 


,(i-l) 


,  ...  ,  d 


0  ’  ”1  ’ 
(i-1) 


That  matrix  of  coefficients  is  a  consistent 


1  '  '  m 

estimate  of  2  times  the  inverse  of  the  covariance  matrix  of  the 
asymptotic  normal  distribution  of  the  as3miptotically  efficient  estimates. 
[See  (1.10).] 

Theorem  3.1.  If  f(X)  given  by  (1.31)  is  positive,  then 

IT 

j.  .  i  r 

"  5  =  —72  I 

p-Xx> 


(3.58)  lim  -  tr  E“^  =  — ^  f  -5^  dX  , 
P  (Ziry  -'-n  f  (A) 
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(3.66)  «  A  5,  ^  Sh  Vh  4  “  i  5,  i  =h 


s,h=.o  -  -s  -  -■>  8  6  -  g^o  -  ~s 


for  any  symmetric  matrices  Gq,  ...  ,  . 

Proof.  Let  F  be  a  nonsingular  matrix  such  that 


8  b 


(3.67) 

(3.68) 


L  =  , 

A  =  FAF'  , 


where  A  is  diagonal  with  diagonal  elements  6--,  ...  ,6  .  [See 

~  11  pp 

Problem  30,  Chapter  6,  T.  W.  Anderson  (1971)  or  Theorem  3,  page  341, 
T.  W.  Anderson  (1958).]  Then 


(3.69) 


Let  H  = 


(3.70) 


0  <  6,  .  <  1  , 

—  ii  — 


i— 1,  ...  ,  p 


m 


z  G  .  Then 
8=0  8~8 


m 


)  tr  L  G  LG,  z  z 
gtiS-O  -  -8  ~  ~h  g  h 


=  tr  FF'  H  FF'  H 


=  tr  F'  H  FF'  H  F 


=  tr  K 

-  t  > 


i>j=l 


where  K  =  F'  H  F  .  Similarly 
(3.71) 


m 


8,b=0 


„  z  z, 

g - h  g  h 


tr  F  A  F'  H  F  A  F'  H 
tr  A  F'  H  F  A  F'  H  F 
tr  A  K  A  K 


-  i!; 


i,j=l 


'ii  "iJ  'ii  • 
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Tf 


1-1  -1  1  f 

(3.59)  lijn  -  tr  S  G  S  =  - » 

P  -  ~g  ~  (27r)2  J--n 


cos  Xg 


A  dX  p  .oo  y  t& 

(2ir)"  ^-ir  r(X) 


(3.60)  lim 

p-XJO 


im  -  tr  G  G  =  —--3-  f 

-g  -  ~h  (2Tr)2  -'-IT 


cos  Xg  cos  Xh 

f^(A) 


dX  , 


g,  h=l,  ...  ,  m  . 

Proof.  The  proof  is  along  the  lines  of  the  proof  of  Theorem  10.2.7 
of  T.  W.  Anderson  (1971).  The  spectral  density  f(X)  is  continuous. 
Therefore,  for  arbitrary  e  >  0  there  exist  autoregressive  processes 
with  covariance  matrices  and  and  (positive)  spectral  densities 

f  (X)  and  f,,(X)  ,  respectively,  such  that 

Li  N  U 


(3.61) 

(3.62) 


fj^(X)  <  f  (A)  <  fy(A)  , 


fL(A) 


fu(A)  - 


<  e 


-IT  <  X  <  IT  , 


-IT  <  X  <  IT  . 


Then,  by  Lenuna  10.2.6  of  T.  W.  Anderson  (1971), 


(3.63) 

(3.64) 


x'  X  _<  x'  Z  X  £  xV  Ey  X  , 


Z7^  X  , 


»./li  '  — 


for  all  X  ;  here  Ej^  and  E^  are  the  p  x  p  covariance  matrices 
corresponding  to  f  (X)  and  f,T(A)  ,  respectively. 

Li  U 

Lemma  3.1.  If  for  A  positive  semidefinite  and  L  positive 
definite 

(3.65)  x'  A  X  <  x’  L  X 

«V  *v, 

for  all  X  ,  then 
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Since  6^^  ,  (3.71)  is  less  than  or  equal  to  (3.70).  This  proves 

the  lemma. 

It  follows  from  the  lemma  that 


(3.72)  ^  tr  G  Y  Yi,  <  5II  tr  G  E~^  G,  y  y, 

irr-^  ~U  ~g  ~U  .  g  •'h  —  ~  ~g  ~  ~h  •'g  •'h 


gth^O  g  ~  ''g 


,  ■ 

<  > _  tr  G  E^  G  y 

-gtl?0  -8--L  -h’'g>'h 


for  every  (yQj  y^»  •••  »  yjj^)  •  Then  for  B  defined  by  (87)  of 
Section  10.2  of  T.  W,  Anderson  (1971)  with  T  replaced  by  p  , 


(3.73)  T~  tr  E"^  G  E“^  G,  y  y,  =  y~  tr  B’B  G  B’B  G,  y  y. 

gth^O  -h^g^h  - g - h^g", 


n  'iz  b  b  g^f  h  .  b  g^^^  y  y 


m,n=l 


fz  (tl  tz 


j,q=l  \  g=0  k,  =1 


tz  fz 

j,q=K+l  \  g=0  k,il=l 


a 


^g)' 


j,q=K+l 


/  ^  K  .  ^  \  2 

(  ^  ®i-k'q-il'  ^g  ) 

g=0  k',£’=0  ^  Jkqjc  gy 


where  K  is  the  order  of  the  approximating  autoregressive  process. 
We  have 


(3.74) 


g^°^  =  1 
®ab  ^  ’ 


a  =  b  , 


=  0  , 


a  5^  b  , 
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(3.75) 


g<h)  =  1 
®ab  ^  ’ 

=  0  , 


a=b+h,  h>0  , 
a?^b+h,  h>0  , 


where,  a,  b,  and  h.  are  Integers.  Then 
8ab  ^  J_^  ®  ^  ’ 

If  5^  =  ^  for  h=0  and  6.  =  1  for  h  >  0  ,  then 
n  2  ti 


/•s  -?Q\  (^)  r  l^(a-b+h)  iX(a-b-h),  ,,  i,  n  i 

(3.78)  ^  2^  I  I-®  '  +  e  '  .  JdX  ,  h=0,  1 . m 


h  >  0  . 


/ 

^-ir 


Then  l/(p-K)  times  the  right-hand  side  of  (3.73)  is 


(3.79) 


( 5Z  e; 

'  ^  j,q=K+l  \  g=0  k’,£*=0 


'’k'  '‘V  h 


6  ,-TT 


^  J"  i  iX[(j-k’)-(q-£0+g]  ^  ^iX[(j-k’)-(q-)l')-gjl  ^5^  y 


2tt  (p-K) 


^  „iX(j-q)  ^ 


j,q=K+l 


xXZ 


iXg  „-iXg 


k’,A’=0 
2 


’’k'  ^ 


-iXk' 


0  1®  0  <e  ®  +  e 

^  8  ■ 


}") 


_ 1 _  (2  ^  y  6  .  .P°s  ^.8 . 


277  fy(X) 


dX 


because 

(3.80) 


"  „_i  T —  K  ,  iXk|2 
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Then  the  limit  of  (3.79)  is 


(3.81) 


4  , 


m 


,1T  ,TT  p-K  . 

I  ZZ„  y,  y,  \  I  I 


(A+V)  (j-q) 


(2tt)"  g,h=0  ®  ^  "  n-x® 


cos  Ag  cos  Vh 

£  /*V  \  £  /-A  - 


f^(A)  f^(v) 


(2tt) 


m  .IT  ^TT 

T  y  y  Yu  ^  I  ^  (^+’^) 

^  gfiS-Q  S  S  1>  ^  P-K 


cos  Ag  cos  vh  ,,  , 

^  nf  ^  A'V  dA  dv  , 


fy(A)  fy(v) 


where 

(3.82) 

Then  (3.81)  is 


k  ir  (^+^)  = 
p-K 


sin^  -I-  (A+V)  (p-K) 
2Tr(p-K)  sin^  j  (A+v) 


(3.83)  ^  y  6  6  11m  ^  k  (A-p)  ^ 

(Zirr  g,h=0  ®  ^  ®  ^  p^  ^-TT  -‘-■n  P  ^ 


m 


- T  /  y  yi.  ^  *5, 

(277)^  g,h=0  g  ^  g 


./: 


cos  Ag  cos  Ah 


dA 


-TT  fj  (A) 


by  substitution  V  =  -y  and  the  argument  leading  to  Theorem  8.3.3  of 
T.  W.  Anderson  (1971) .  Then 


m 


(3.84) 


lim  1  V —  r>“l  T.— 1 

— —  )  y  y.  tr  E  G  E  G, 

P  P  g'^o  ®  ^  ~  ~g  ~  ~h 


is  greater  than  or  equal  to  (3.83).  Similarly 
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(3.85) 


m. 


lim  ^  LI  Yg  Yh  tr  G  G 
p-x»  P  g,h=0  ~  ~g  ~ 


,  m  ^ir 

-  ^  g5o  L 


cos  Ag  cos  Xh 

4 '« 


dX 


Since  e  is  arbitrary 


m 


(3.86)  lim  —  }  y  y,  tr  Z  ^  G  Z  ^  G, 

p-^  P  gTh^o  S  h  .  .g  .  .h 


m 


Yi.  <S„  6. 


(27r)^  g,h=0  g  S 


IT 

h,  f 

J-TT 


COS  Xg  cos  Xh 
(X) 


Since  (3.86)  holds  for  every  set  of  yQ,  y^,  ...  *  Y^^  >  the  theorem 
follows . 

It  is  of  interest  to  compare  this  result  with  Theorem  8.3.3  of 
T.  W.  Anderson  (1971)  which  gives  the  asymptotic  covariances  of  the 
sample  covariances  of  a  stationary  process  with  a  continuous  spectral 
density.  The  spectral  density  of  an  autoregressive  process  with  co¬ 
efficients  a^,  . . ,  ,  is  proportional  to  the  reciprocal  of  f(X) 
for  the  moving  average  process;  hence  the  Gaussian  part  of  (37)  of 
Section  8.3  of  T.  W.  Anderson  (1971)  is  proportional  to  the  coefficient 

of  y  y  in  (3.86). 

§  ^ 
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4.  Asymptotic  Theory  as  the  Sample  Size  Increases 

When  . Xj^  are  observed  with  ^x  =  y  satisfying  (1.2) 

and  Z  is  knovm,  the  Markov  estimates  of  3^,  ...  ,  are  the 
solutions  to 


(4.1) 


i=l 


j=l,  ... 


In  matrix  notation  the  solution  is 


(4.2) 


!n  -  '5'  r'  r  r"  > 


where  “  (^2’  *“  ’  ^r^  ' 


(4.3) 


Z  =  (z. ,  ...  ,  z  )  . 
~  ~r 


(We  use  the  subscript  N  on  x^^  and  8^^  to  emphasize  the  dependence 
on  N  . )  Then  the  covariance  of  (3jj-*8)  is 


-1  „^-l 


(4.4)  P[/N  (3jj-3)J  “  N  ^(3jj-3)(3jj-3)’  =  (Z’  Z) 


Regardless  of  the  distribution  of  X  ,  i^N  (3jj-3)  has  a  limiting  normal 

distribution  with  mean  vector  0  and  covariance  matrix  (4.4),  because 

/tf  (Xjj-y)  has  a  limiting  normal  distribution. 

If  Z  is  unknown,  let  us  suppose  that  we  have  a  consistent  esti- 

mate  Z„  of  Z  .  Consider  the  estimate 
~N 


(4.5) 


is  =  '5'  S"  ?>■'  r  S'  iK  • 


Then 
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(4.6) 


=  [(Z'  E"^  Z)~^  Z»  -  (Z'  E"^  Z)"^  Z*  E“^]  jM  (Xjj-Z3) 


converges  stochastically  to  0  because 


(4.7) 


plim  (Z*  E"^  Z)"^  Z‘  ij-  =  (Z’  E"^  Z)"^  Z’  E"^ 


N 


and  /N  (Xj^-Z3)  =  ^  limiting  distribution.  Thus  /N  (3^-3) 

has  a  limiting  normal  distribution  with  mean  0  and  covariance  matrix 
(Z'  E  ^  Z)  ^  ,  which  is  the  same  as  the  limiting  normal  distribution  of 

A  A 

3j^  .  If  3jjj  is  asymptotically  efficient, then  3^  is  (in  the  same  sense). 

A 

In  particular,  when  3j;j  is  maximum  likelihood,  as  when  X  has  a  normal 
distribution,  it  is  asymptotically  efficient  in  the  sense  of  attaining 
the  Cramer-Rao  lower  bound  for  the  covariance  matrix  of  unbiased  estimates. 
We  summarize  the  result. 


Theorem  4.1.  Let  x^ ,  ...  >  be  Identically  distributed  with  mean 

A 

^  X  =  Z3  and  covariance  matrix  E  and  let  E„  be  a  consistent  estimate 

A  A 

of  E  .  Then,  if  3^  is  given  by  (4.5),  /n  (3^t“3)  has  a  limiting 

^  ■  ~jM  's/JN 


normal  distribution  with  covariance  matrix  (Z  E  ^  Z)  ^ 

A 

by  (4.2)  is  asymptotically  efficient,  so  is  3^  . 


If  ^y.en 


We  now  apply  the  result  to  the  estimation  of  0^,  O^,  ...  >  • 

Theorem  4.2.  Let  x^,  ...  >  x^^j  ^  N  observations  from  N(ii,  E)  , 

where  y  is  known  and  E  is  given  by  (1.3).  Let  C  be  defined  by 

(1.7).  Let  ...  >  be  a  consistent  set  of  estimates 

-  .u  1  m  • 


of^  Oq,  ...  ,  0^ 


Let  ...  ,  be  the  solution 

-  0  1  m  - 
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tg  (2.1)  for  i=l  .  Then  (a.^^^-cf,),  ...  ,  (d^^^-O  ) 

-  -  u  u  11  m  m 

have  a  limiting  normal  distribution  with  means  0  and  covariance  matrix 
(1.10),  and  •••  >  are  asymptotically  efficient. 


Proof.  In  Theorem  1  replace  by  C  ,  Z3  by  )  > 

S  by  $  ,  and  Z  by  inserting  •••  >  in  $  . 

Then  is  a  consistent  estimate  of  $  and  the  limiting  normal 

distribution  follows.  Since  the  solutions  to  (1.22)  are  as3nnptotically 
efficient,  the  theorem  follows. 

If  u  is  unknown,  it  can  be  estimated  by  ,  and  C  in  Theorem  2 

can  be.  replaced  by  (1/N)  ) _ _  -  (x^-x„)  (x^-x.,) '  .  If  y  has  the  form 

- 01“1  ~ct  ~ 

(1.2),  it  can  be  estimated  by  )  j  i  3.  z.  ,  where  3,,  ...  ,3  is  a 

^ - j=l  j  ~3  ’  1’  ’  r 

solution  to  (2.21).  For  the  asymptotic  theory  uses  the  fact  that 
.  (C-E)  has  a  limiting  normal  distribution, 

1  N 

(4.8)  [C  -  -  YZ  ’ 

a=l 


(4.9)  A  [C-C]  =  /n  (x^-u)  (^-uy  -  »¥  (Xj^-C)(Xj^-a)’ 


and  the  facts  that  plinij^  ^  =  y  ,  plinij^  ^  =  y  ,  and  y¥  (x^^j-y) 
and  /n  have  limiting  (normal)  distributions. 
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5.  Estimation  of  the  Coefficients  of  a  Moving  Average  Process 

We  now  treat  the  problem  of  estimating  the  coefficients  a^,  a^, 

2 

of  the  moving  average  process  (1.25)  with  the  restriction  a  -  1 


replacing  the  restriction  aQ  =  1.  Then 

m-h 

h  ^  g  g+h 


(5.1) 


h=0 ,  1 ,  ...  ,  m 


Let 

(5.2) 

Then 

(5.3) 


m 


M(z)  =  )  a.  z 

j=o  J 


m-j 


m 


^  =  M(z)  M(z"^) 

h=-m 


by  (5.1). 

Let  the  roots  of  M(z)  =  0,  which  is  (1.32),  be  z, ,  . . .  ,  z  . 

1  m 

(5.4) 


m 


M(z)  =  a-  YJ  (z-z  ) 
j=l  ^ 


and 


m 


/  c  c  \  \ ^  m+h  2  m 

(5.5)  2 _  CTjj  z  -  ^ 

h=-m 


m 


TT  (z-z.) 


j=l 


a. 


m 


TT  (z-z . ) 

,  j=l  J 


1 


r 


m 


TT  (l-s;  z,) 

L  j-1  ^ 


Thus,  if  a(m)  0,  the  2m  roots  of 
(5.6) 


m 


h=-m 


^  m+h  « 
Oh  Z  -  0 


2-  j  «  «  *  )  9  1/z  • 

•L  in  X  m 


Then 
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Conversely,  if  Oq,  0^,  ...  ,  0^  0  are  given,  the  set  of  2m 

roots  of  (5.6)  can  be  partitioned  into  two  sets  of  m  roots  each,  such 
that  the  roots  in  one  set  are  less  than  or  equal  to  1  in  absolute  value 
and  are  the  roots  of  an  equation  (1.32)  with  real  coefficients,  that  is, 
M(z)  =  0,  defined  by  (5.2)  with  real  coefficients,  and  the  roots  in  the 
other  set  are  reciprocals  of  the  respective  roots  in  the  first  set.  The 
roots  define  the  coefficients  in  M(z)  except  for  normalization,  which 
is  determined  by  (5.1)  for  k  =  0.  [See  Section  5.7  of  T.  W.  Anderson 


(1971).] 

The  equation  (5.6)  can  be  written 

m 


(5.7) 


A  \ ^  ,  m-hv 

h=l 


m 


=  z 


m 


2b  ,  , 

h»l  z 


Let 

(5.8) 


w  =  z  + 


1  z  +  1 


that  is, 

(5.9) 


zw  =  z  +1 


Then 


TAX  2r  2r  .  ,.2r  ^  2  x2r 

(5.10)  z  w  =  (zw)  =  (z  +1) 


=  (z2)2’'+ll 


+  2r  z^  (z^)^’^"^  +  1 


•  •  •  +  ('ll)  [(=■')'  +  l]  +  f 


2r 


A3 


2r+l  2r+l  ,  .2r+l  ,  2.,.2r+l 

(5.11)  z  w  =  (zw)  =  (z  +1) 


n  2.2r+l  ^ 
=  (z  )  + 


^  +  (2r+l)  ^ 

H.  l]  +  ...  +  (z^)'  p«] 


2,s 


[The  coefficients  of  (z  )  +1  in  (5.10)  and  (5.11)  are  the  same  as 

the  coefficients  of  A  in  (16)  and  (17)  of  Section  6.5  of  T.  W.  Anderson 

s 

2  s 

(1971).]  When  we  solve  (5.10)  and  (5.11)  successively  for  (z  )  +1,  s=l, 

2,  ...  ,  m,  we  obtain 


2,2r, 


2^2r+l 


2r 

2r  . 

2r-2 

•  - 

=  z 

[“  ‘'2r,2 

w  +  . . .  +  c 

2r,2r 

9 

2r+l  r  2r+l 

2r-l 

+  i  =  z  ^w 

+  ^2r+l,2  ” 

+  » •  < 

'  ‘^2r+l,2r 

’  ^2r,2r’  ‘^2r+l,2’ 

•••  ’  '^2r+l,2r 

are  appropriate 

w 


constants.  Then  (5.7)  can  be  written 


(5.14)  0  =  {Oq  +  O^w  +  O2  [w^-2] 


2n 


■  ■  r  I  2xx**‘*2  I 

+  ...  +  a2n  [w  +  C2^^2  ^  +  ...  +  C2j^^2n 


]} 


=  z  "  {02^  w  +  w  ^  +  ...  +  Qq  -  202  +  ...  +  02^  ^2n,2n^ 


m=2n 


(5.15)  0  =  {Oq  +  a^w  +  02  [w^-2] 

.  ,  -  r  2n+l  2n-l  ,  i 

+  •••  +  “2n+l  +  %+l,2  ”  +  •••  +  <=2,ri-l,2n  > 

2n+l  2n+l  ,  -  2n  ,  ,  ^ 

"  ==  ^'^2n+l  +  '^2n  +  *^0  ~  ^'^2  ''2n  ''2n,2n^  ’ 


m  =  2n+l  . 
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Solving  (5.14)  or  (5.15)  as  an  equation  in  w,  we  obtain  m  roots. 


w, ,  . . .  ,  w  .  Then  solve 
1  ’  m 


(5.16) 

for  z,  to  obtain  the  roots 

(5.17) 


z  -  w^  z  +  1  =  0 


w. 

_i  +  /(^)2  _ 
2  -  \2  I 


to  (5.7);  the  pair  of  roots  to  (5.17)  are  reciprocals.  The  m  roots 

with  absolute  value  less  than  or  equal  to  1  are  the  desired  roots  of 

M(z)  =  0  and  the  coefficients  of  the  polynomial  M(z)  are  the  desired 

ttQ,  a^,  ...  ,  except  for  a  constant  of  proportionality. 

We  want  to  modify  the  numerical  procedures  discussed  in  Section  2. 

If  ...  ,  are  the  estimates  at  the  i-th  stage, 

let  ...  ,  be  the  desired  solution  of  (5.1)  with 

0  1m 

the  replaced  by  the  ^^'s.  We  want  to  determine  the  next 


iteration;  let 


(5.18) 


a. 


(i) 


(i-1) 


+  d 


h’ 


h=0,  1,  ...  ,  m 


Then  if  we  substitute  into  (2.1)  with  C  =  xx*  we  obtain 

m-f 


m 


-1  .  ;-l  .  ^  ,.(i-l)  ^  ,  /^(i-1) 


(5.19)  XZ  tr  G  lit,  G.  YZ 

Y-  n  -L  ~8  -L  ~t  V  TT 


f=0 


h=0 


''-1  ''-I 

=  x’  E,  T  G  E/,  X  , 
~  ~i-l  ~g  ~i-l  ~ 


‘^n^  Kir  +  Vf> 


g—  0,1,  ...  ,m. 


The  left-hand  side  of  (5.19)  can  be  written 
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(5.20)  ^  tr  5-1,  0  5-1,  ^  +  a<if>  d^  +  d^^  d^)l 


=  tr 


1-1  ~g  ~1-1  ~f  -  f 

G  +  iz  tr  G  i:\  g, 

tr  E“^  G  ET\  G  d, 

^  ^i-1  ^g  ^1-1  ..k  j+k  J  : 


ni  ^1  ^1  in— £ 

)  tr  ZT  G  ZT  ,  G.  y~  d,  d,  _  . 

'•  ■  ■■  ~i-l  ~e  ~i-l  ~f  t — 7T  h  h+f 
f=0  h=0 


Using  the  linear  terms  only,  we  can  obtain  from  (5.19)  the  equations 


(5.21) 


m  f- 

n 

j=0  L-k=0 


^  „  ;-l  ^  .^(i-1) 

)  tr  2,,  -  G  Z.  T  G,  a;  , 


m  .j  /s—.  /\.| 

+  g  tr  G^  s'-  > 

=  x'  Z  G  Z  .  X  -  tr  Z.  _  G  , 

-  ~i-l  ~g  ~i-l  ~  ~i-l  ~g 


g— 0 ,  1,  ...  ,in. 


If  •••  j  are  consistent  estimates  of  a™,  a,,  ...  ,  O  , 

0  1m  0  1  ’  m’ 

then  ...  ,  are  consistent  estimates  of  a^,  o,. . a  . 

O’l’  ’m  0’l’’m 


If  the  solution  to  (5.21)  for  1  =  0  is  d^,  d,  ,  ...  ,  d  ,  then  a^^^  +  d^. 

0  1m  0  0 

af^^  +  d- ,  ...  ,  a^^^  +  d  are  consistent  estimates  of  a„,  a,,  ...  ,  a  . 
11mm  .  0  1  m 

It  may  be  expected  that  these  estimates  are  asymptotically  efficient  as 

p  00  . 

A  method  suggested  by  Durbin  (1959)  for  estimating  a^,  ...  »  ci^ 


when  ttQ  =  1  is  to  solve 


(5.22) 

for  some  n  >  m. ,  where 


n-1 

'>  c.  ,  -9.  =  -  c. 

J=o  J 


•L^X  «  «  •  •  a  XI 
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(5.23) 

1  ^ 

j  -j  P  ^  1  i+3 

j*=0,  •••  jXij 

and  let 

=  1.  Then  solve 

(5.24) 

m 

)  h,  6i  =  -h,  , 

f-S  s  i 

}  m  J 

where 

(5.25) 

*'g  '“-s'  ^t+g  ■ 

9  9  m  • 

See  T.  W.  Anderson  (1971),  Section  5.7.2,  for  more  explanation  and  dis¬ 
cussion.  Raul  Mentz  has  shown  that  the  resulting  estimates  are  approximately 
consistent  and  asymptotically  noraally  distributed  as  p  These  esti¬ 

mates  (suitably  normalized)  could  be  xised  as  starting  values  in  the  above 
iteration  procedure. 

Since 

(5.26)  plim  c  =  f  ,  j=0,  1,  ...  ,  m  , 

p-XXJ  J  J 

(5.27)  plim  c,  =  0  ,  j=m+l,  ...  ,  n  , 

p-xx)  ^ 

the  equations  (5.22)  can  be  replaced  by 

(5.28)  C*  t  =  c  , 

where  C*  has  the  form  of  (1.29)  with  O  replaced  by  c  ,  g=0,  1,  ...  , 

S  8 

nij  ?  =  (?-j^»  •••  »  ^  •••  »  forward  solution 

can  be  applied  to  (5.28)  and  the  forms  h.,  ...  ,  h  can  be  formed  as 

U  m 

indicated  in  Section  2  with  p  replaced  by  n. 
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